home *** CD-ROM | disk | FTP | other *** search
-
-
- { -> 278 }
- procedure simps(function f(x: real): real;
- lower,upper,tol : real;
- var sum : real);
-
- { numerical integration by Simpson's rule }
- { function is f (as parameter), limits are lower and upper }
- { with number of regions equal to pieces }
- { partition is delta_x, answer is sum }
-
- var i : integer;
- x,delta_x,even_sum,
- odd_sum,end_sum,
- end_cor,sum1 : real;
- pieces : integer;
-
- function dfx(x:real):real;
- begin
- dfx:=-1.0/sqr(x)
- end;
-
- begin
- pieces:=2;
- delta_x:=(upper-lower)/pieces;
- odd_sum:=f(lower+delta_x);
- even_sum:=0.0;
- end_sum:=f(lower)+f(upper);
- end_cor:=dfx(lower)-dfx(upper);
- sum:=(end_sum+4.0*odd_sum)*delta_x/3.0;
- repeat
- pieces:=pieces*2;
- sum1:=sum;
- delta_x:=(upper-lower)/pieces;
- even_sum:=even_sum+odd_sum;
- odd_sum:=0.0;
- for i:=1 to pieces div 2 do
- begin
- x:=lower+delta_x*(2.0*i-1.0);
- odd_sum:=odd_sum+f(x)
- end;
- sum:=(7.0*end_sum+14.0*even_sum+16.00*odd_sum
- +end_cor*delta_x)*delta_x/15.0;
- until abs(sum-sum1)<=abs(tol*sum)
- end; { simps }
-
-
-